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THE BALLISTIC PEOBLEM 
By Frank Gilman 

The problem of the motion of a projectile in a medium resisting as any 
power of the velocity was first solved by John Bernouilli in 1718, but he ex- 
pressed the coordinates of any point of the trajectory in a form where a final 
integration is required to determine them numerically ; and, as far as the gen- 
eral problem is concerned, no later mathematician has improved on his results. 

Modern experiments in gunnery, notably those conducted for the British 
government (1865-1880) by General Francis Bashforth, have established the 
fact that for modei-ate velocities, i. e. from 100 to 900 feet a second, the re- 
sistance of the air to the passage of the projectile is very nearly proportional 
to the square of the velocity ; that for greater velocities i. e. from 900 to 1300 
feet a second it is nearly proportional to the cube of the velocity ; and that for 
high velocities, i. e. from 1300 to 2780 feet a second it is again very closely 
proportional to the square of the velocity. 

In 1873 General Bashforth published a report of his experiments accom- 
panied with the mathematical theory for the case where the resistance of the 
air is as the cube of the velocity, and an elaborate set of tables based on the 
method of quadratures for the integration of the final equations. 

In 1890 he published a revised report in which the more difficult case 
where the resistance is as the square of the velocity was treated in like manner. 

The same subject has been somewhat differently treated by Niven,* 
Siacci,t and Ingalls.t 

In this paper I shall attempt to show that by means of the application of 
the method of Least Squares to the development of functions the diflFerential 
equations for the trajectory when the resistance is as the square of the velocity 
can be integrated without excessive labor, so that they will give results of 
very considerable exactness. 

I shall first consider a method oi very great accuracy but requiring much 

* Niven, Proceedings of the Royal Society, London, Vol. 26, 1877. 

t Siacci, Balistique Exlerieure; Traduction annot^e par P. Laurent; Paris, 1892. 

% Ingalls, Exterior Ballistics, Washington, 1900. 
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labor in its application, and I shall afterwards develop much simpler working 
methods of sufficient accuracy for practical purposes. 

The following* are the well known differential equations of motion for 
the case under consideration : 
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(2) 
(3) 



cdx 




dp 








pyjl + p^ 


+ logCiJ + 


V/1 + p-') - 


-0" 


cdy 




pdp 






Psli+P"" 


+ log(i> + 


v/i+y)- 


■ c' 


cn.dt. 




— dp 







^0-p\ll-\-p^-\o^{p+ VI +y) 



in which dx denotes the differential of the abscissa and dy the differential of 
the ordinate, while ^, or -p, is the tangent of the angle at which the trajec- 
tory is inclined to the horizontal. 
We have also 



^ = 2cAcos^0o "^ ^"^^ +^" + ''*^^^" + ^^ +^"^ 



in which c=j^,g being the acceleration of gravity and k the velocity of the 



v\ 



projectile when the resistance becomes equal to its weight ; h — -^ where Vq 

is the initial velocity ; <^o is the initial angle of elevation and p^ the initial value 
of^. It is easily seen that O is positive and always greater than 

p\lTTp^ + log(iJ + sliTp'^.) 

First metliod. Let us replace the transcendental function p^\ + p^ 
+ log(^ + v^l + /J*) which we will represent byy(^) by the rational expression 
a'p + b'p^ + c'p^ in which the coefficients a', U, and c', are to be calculated by 
the method of Least Squares, the general form of the equation of condition be- 
ing as follows : 

a'p + b'p^ + c'p^-f{p) = 0. 

Make successively in this equation^ =; 0. 1,^ = 0.2, etc. upto^ = 1.0, 
* V. Bartlett's Analytical Mechanics. 
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thus forming 10 equations of condition. Forming the normal equations and 
solving them we obtain 

Log a' =0.3012124 

Log 6' = 9.5126277 
Log (-c') = 8.4899010 

A more detailed exposition of this method of developing irrational and 
transcendental expressions is given in an article that was published in the 
Journal of The FranMin Institute for June, 1894, entitled "Application of the 
Method of Least Squares to the Development of Functions." 

The above values of a', h' , and c' when substituted in the rational expres- 
sion will make it give values for_/'(j9) whose greatest eiTor will hardly exceed 
one in the fourth place of decimals, while the integrals, ex and cy, are likely 
to be even moi-e exact. 

In formulas developed by the method of Least Squares (especially when the 
number of conditional equations is large) the sum of the positive errors in a 




series of results, will be nearly equal to the sum of the negative errors 
Therefore if the line Oafd represents the course of a given curve with reference 
to the axes x and y, its location as determined by the method of Least Squares 
is shown by the curve Ocd ; and although these curves may deviate sensibly 
from each other at certain points, as at a and c, yet the areas included between 
the lines Oe, ed and each curve will be practically the same. In other words 
the corresponding definite integrals will be nearly equal. 

The differential equations (1) and (2) may now be written as follows : 

^^ (4) 



cc'dx = 



h' a' C" 

-^ e ^ d d 



eddy = 



pdp 



^ b' a' C 



(5) 
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These equations can now be integrated by the method of decomposition 
into partial fractions. 

h' a' O 
Let I = -„ ill = -,, n = J and observe that I and m are negative and 

CO c 

n positive. 

By Descartes' Rule of Signs p^ + Ip^ + mp + n = has two positive real 
roots, one negative real root and two imaginary roots ; call them a, fi, 7, 

a 4_ e v/~l, S - e ^~l. Then 

1 _ L M IT Pp+Q 

+ z. — a + ;; — z + 7Z — SWTT2 ' (*') 



p^ -\- Ip^ ->r mp -{■ n p — a p — l3 p — y (P — ^Y^ + 

L = r p~"- 1 = ^ 

\_p6 + lp^ + mp + n]p^^ 5a* + ila^ + 

1 



m 
Similarly, 

M = 



N = 



50* + Up^ + m ' 

1 
5/ + 3;,y2 + m ' 



Clearing ((5) of fractions and writing the coefficients of^/ andjr;^ equal to zero, 
noting that a + /3 + 7 + 28 = 0, we have 

L + M + N -^ P = i), &nA La + M0 + N'y + iPh Jr Q = 0, 

whence 

P=-(Z + il/ + iV;, 

Q = -(La + i»f/8 + JV7 + 2P8). 
E(juations (4) and (5) now integrate into 
vc'x = L\og{p - a) + ili log(^ -/3) + J^\og(p - 7) + ^ \ogl{p - B)^ + e'^'] 

cc'i/ = La log {p- a) + M0 log (i? - /3) + JSTy log {p - 7) 

The constants of integration are determined from the condition that when 
p is ecjual to its initial value p„, cb = and t/—0; but the formulas may be used 
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as they stand without the determination of the constants of integration if p is 
taken between proper limits. 

We will now apply these equations in solving the following numerical ex- 
ample which is given by Professor Bashforth in his treatise on the Chrono- 
graph, page 97. 

Initial angle of elevation <^o = 30° ; initial velocity F= 751 ft. per sec. ; 
diameter of shot = 6.27 in. ; weight = 701b3. Here c = 0.0000339822, 
(7=3.455425. 

We first find the five roots of the equation 

p'^ + Ip^ ■{■ mp + n = 0, 

in vvhichZ = ^=- 10.5372354, m=^ = -64.7606866, n = --?= 4-111. 84085. 
c' c d 

The following are the values of the roots : 
a = + 1.3779718 
/3=: + 3.6221343 
7 = - 4.0293468 

S zz= - "" ^f + '^ = - 0.4853796 

£« = 38* + Z - 7(a + iS) - a/3 = + 5.3255100 

Next computing the five constants by the formulas of the preceding page 
we get : 

X = - 0.0093670, il'/:= + 0.0026236, i\r= + 0.0013514, 
P = + 0.0053920, § = + 0.0140848. 

It is evident that the accuracy of the integrals as given by equations 
(7) and (8) will depend on the exactness with which the following equation 
is satisfied : 

1 



PV^l ^f + \o%{p + v^l +y) - C 

e\j) — a p — p P — "i (i> — o) + e^J 
Substituting in this e(|uation the values of the roots and constants already 
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found, wo obtain the following tabic, which gives the values of the two mem- 
bers of this equation corresponding to ditl'ercnt values of j? : 

p naar membek 

O.Oo 0.2!»803 

0.1(» 0.30721 

0.15 0.31703 

0.20 0.327.'>7 

0.2') 0.33805 

0.30 0.35130 

0.35 0.36478 

0.40 0.37957 

0.45 0.39590 

O.-'iO 0.41404 

0.55 0.43436 

0.60 0.45727 

0.65 0.48334 

0.70 0,51327 

0.75 0,54801 



SECOND MEMBER 
0.29803 

0.30722 
0.31704 
0.32758 
0,33896 
0.35132 
0,36479 
0,37957 
0.39589 
0,41403 
0.43434 
0.45725 
0.48331 
0,51325 
0,54801 



In accordance with what has been said in connection with figure 1, it 
follows that the integrals are likely to be even more exact than might be in- 
ferred from the indications of the above table. 

Substituting the values of the roots and constants in equation (7), and 
taking it between the limits ^ = tan 30° and p = 0, we find x = 6026.8 feet, 
which is the horizontal distance from the origin to where the projectile attains 
its greatest elevation. 

Professor Bashforth finds by (quadratures x = 6032 ft,, which is 5 feet in 
error. 

Taking equation (8) between the same limits, we obtain y = 1864.4 ft. 
Professor Bashforth finds y = 1867 ft. In order to find the range we will as- 
sume with Bashforth that the angle of fall = - (37" 3'), and taking ecjuations 
(7) and (8) between the limits jj = tan 30°and^ = tan — (37° 3'), we obtain 

X = 11326.1 ft. = range, 
y = - 7.2 ft. 

According to Bashforth, x = 11334 ft., and y = — 6,5 ft., showing an error 
f>f about 8 foot in the ranjre. 
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The equations (7) and (9) undoubtedly give the coordinates of* the tra- 
jectory correct to the nearest foot ; but their use involves great labor, and 
they are chieily valuable as a criterion by which to compare the results of 
simpler methods. 

Second method. If instead of replacing our f{p) =p^l + p^ 
+ log(i> + VI +i>') by a'p + ftjp* + dji^ we assimie f{p) = op + 6p*, a and 
b being constants to be determined as before, our method while open to the 
objection that two separate sets of formulas are required one for the ascending 
and the other for the descending portion of the trajectory, is very much less 
laborious in practice than the one we have been considering. 

Since /(—i>) = —f(p), if we determine a and b on the assumption that 
p is positive we must use ap — bp^ iovf{p) when j? is negative. 

Integrating the diflTerential equations (1), (2), and (3) on the assump- 
tion that/( j>) = op + 6p* we get for the ascending portion of the trajectory : 
i , g + 2/^ - yg" + ^bV 



X = 



cy/a* + ibC ^ a + 2bp + yja* + ibC ' 



y=^Aog{C-ap-b]^)-ff. (11) 



1 . , / g + 26i> X 
For the descending portion where /(^) = ap — bp* Ave get : 

X = ^ log « - 26p - y/^-^ - ^^'t.' ^ .^^, ^ 4^^,^ 

cy/a^ -ibO a- 2bp + ^a^ - ibO 

X = _ .^^ tan-i ( "~^^^ ^ , if a^ < 46C, (13') 

cyJAbV - a* \y/4<»C' - aV 

'■^ = -Wc ^'^SiO-ap + bp") + H , (14) 

t = ^ log[s/0-ap + bp^ - p^b + 2^J . (15) 

As the constants of integration have not been introduced in the above 
tbrmulas, we must remember in using them to take p between the limits for 
which X, y and t are to be calculated. 
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The following 


table gives the values of a and b corresponding to difterent 


igles 


of elevation 










V 


4 


> 


logo 


log h 


a 


a* 


0.1 


b" 


43' 


0.3007161 


8.6654758 


1.998555 


3.994222 


0.2 


11 


19 


0.2998279 


8.9559684 


1.994472 


3.977917 


0.3 


16 


42 


0.2983910 


9.1254675 


1.987884 


3.951682 


0.4 


21 


48 


0.2964379 


9.2440980 


1.978964 


3.916298 


0.5 


26 


34 


0.2940089 


9.3342507 


1.967927 


3.872735 


0.6 


30 


58 


0.2911508 


9.4060248 


1.955018 


3.822097 


0.7 


35 





0.2879118 


9.4649223 


1.940492 


3.765508 


0.8 


38 


40 


0.2843400 


9.5i.n827 


1.924598 


3.704077 


0.9 


41 


59 


0.2804821 


9.5563036 


1.907577 


3.638850 


1.0 


45 





0.2763813 


9.5925171 


1.889650 


3.570776 



In using this table no interpolations will be necessary, as either one of two 
consecutive sets of values of the a and b functions can be used without making 
any material difference in the results. The quantities a and b in the table are 
determined in the usual manner by the method of Least Squares, the general 
form of the equation of condition being 

ap + 6p* -f{p) = 0. 
Make successively in this equation j> = 0.01,^ = 0.02, etc., up top = 0.10, 
thus forming 10 equations of condition. Solving the resulting normal equations, 
we obtain the log a and log b given in the first horizontal line of the table. 
In the same manner the log a and log b found in the second horizontal line of 
the table are fqrmed from 20 equations of condition, and those in the last-hor- 
izontal line are formed from 100 equations of condition. 

We will now solve by the above equations and table the numerical exam- 
ple thatwe have already worked by the first method, and in which O = 3.455425, 
c = 0.0000339822. Since in this example the initial angle of elevation equals 
30°, we take the values of a and b from that line in the table corresponding to 
p = 0.6. Substituting these quantities together with the values of O and c 
in equations (10) (11) and (12) and taking them between the limits corre- 
sponding to ^ = tan 30° and p = 0, we obtain 

x = 6026.7 ft., 

y = 1864.8 ft., 

t = 10.38 sec. 
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As the angle of fall is 37° 3', we take for the descending branch the values 
of a and b corresponding to^ = 0.8. Substituting these in the equations (13') 
(14) and (15) and taking them between the limits j) = tan (37° 3') and^ = 0, 
we find 

X = 5300.7 ft., 

y = - 1871.1 ft., 

t= 11.13 sec. 

Eange = 6026.7 + 5300.7 = 11327.4 ft. 

Time of flight = 10. 38 + 11.13 = 21.51 sec. 

If we use the values of a and b corresponding to p = 0.8 for the ascending 

branch of the trajectory, we obtain the following for the x, y and t of this 

branch : 

X = 6026. 2 ft., 

y= 1865. 3 ft., 

t = 10. 38 sec. 

From this it appears that we may generally use the same values of a and 
b for the ascending as for the descending branch, those quantities being taken 
Avhich correspond to the angle of fall. 

We will illustrate by one more example taken from Capt. Ingalls' Treatise 
on Exterior Ballistics, page 104. 

A 12-inch service projectile is fired at an angle of departure of 10°, with 
an initial velocity of 1886 ft. per sec. ; it is required to find Xi, y,, <j, corre- 
sponding to p — 0, and x^, y-i, t^ for the point where the angle of fall is — 13°, 
both sets of coordinates being referred to the mouth of the gun. 

In this example O = 0.72037, c = 0.0000254985. Take a and b from 
that line of the table corresponding j) = 0. 2. Substituting in equations (10), 
(11) and (12), and taking j9 between the limits tan 10° and 0, we find 

xj = 13203.8 ft., 
yi z= 1294.4 ft., 
ti =8.46 sec. 

For the descending branch take the a and b corresponding toj? = 0.3, 
and substitute in equations (13) (14) and (^15), takings between the limits 
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tan (— 13°) and 0. We thus obtain the following values of the coordinates re- 
ferred to the vertex of the trajectory ; 

x= 9705.3 ft., 
y = -1027.3 ft., 
t = 8-32 sec. 
Whence 

x.i = 13203.8 + 9705.3 = 22909.1 ft., 
y.^= 1294.4-1027.3= 267.1ft., 
t.2 = 8.46 + 8.32 = 16.78 sec. 

Capt. Ingalls working by "Niven's Method" obtains the following values 

for these quantities : 

xi = 13180.7 ft., 

yi= 1291.8 ft., 

a;2 = 23007.0 ft., 

2/2= 241.7 ft., 

<i = 8.44 sec, 

t.i = 16.83 sec. 

If the actual law of the resistance is that of the square of the velocity these re- 
sults of Capt. Ingalls are somewhat erroneous, the x.^ being 98 feet in error. 
The more exact equations (7) and (8) give the following values for these 

coordinates : 

Xi = 13205.6 ft., 
yi= 1295.5 ft., 
x.i = 22908.9 ft., 
y2= 267.3 ft. 

Third method. A rougher method which in point of simplicity leaves 
little to be desired, and yet is surprisingly accurate in application is based on 
the assumption that f(p) =p^l + p^ + log(j3 + y/1 + p'^) = fip. 
E(iuations (1), (2) and (3) integrate at once into 

cx==-log(C-ap). (1(5) 

''^=1 [p-^^^''SiO-ap)'j, (17) 

2 

\Jcff . ( = -y/0-ap, (1«) 
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and the ascending and descending portions of the trajectory do not need to be 
dealt with separately. 

Here again as the constants of integration are not written we must in 
using the equations take the diiferencc between the values corresponding to 
the extreme values of^. 

The value of a is to be determined by the method of Least Squares, our 
equation of condition being ap —f{p) = ; but as we have but one unknown 

•f{P) 



we have merely to compute the arithmetical mean of the values of" 



P 



The following table gives the values of a for different angles of elevation : 



p 




4> 


a 


P 


«/> 


a 


0.05 


2° 


52' 


2.00028 


0.55 


28 


49 


2.03364 


0.10 


5 


43 


2.00124 


0.60 


30 


58 


2.03976 


0.15 


8 


32 


2.00272 


0.65 


33 


1 


2.04635 


0.20 


11 


19 


2.00474 


0.70 


35 





2.05339 


0.25 


14 


2 


2.00731 


0.75 


36 


52 


2.06087 


0.30 


16 


42 


2.01041 


0.80 


38 


40 


2.06879 


0.35 


19 


17 


2.01403 


0.85 


40 


22 


2.07711 


0.40 


21 


48 


2.01818 


0.90 


41 


59 


2.08585 


0.45 


24 


14 


2.02284 


0.95 


43 


32 


2.09497 


0.50 


26 


34 


2.02799 


1.00 


45 





2.10447 



Applying this third method to the examples already worked by the first 
method and the second method, we get for the range and time in the first ex- 
ample 

x= 11335 ft., 

t = 21.51 sec. ; 

and for the coordinates and the time in the second example 

x.i=z 22906 ft., 
y.2= 266 ft., 
ti = 16.78 sec. 
CHURcnvtixR, N. Y. 



